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Abstract 

Modes of speciation have been the subject of a century's debate. Traditionally, 
most speciations are believed to be caused by spatial separation of populations (al- 
lopatry). Recent observations El CO and models 001311113111311111121113 
El El) show that speciation can also take place in sympatry. We discuss a compre- 
. hensive model of coupled differentiation in phenotype, mating, and space, showing 

that spatial segregation can be an induced process following a sympatric differen- 
tiation. This is found to be a generic mechanism of adaptation to heterogeneous 
environments, for which we propose the term diapatric speciation |16j . It explains 
the ubiquitous spatial patching of newly formed species, despite their sympatric 
origin [21|31|I|. 
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Allopatric speciation occurs in populations extending over a sufficient range in space 
and time. If subpopulations become spatially isolated, they can diverge in phenotype by 
adaptation to different environments as well as by genetic drift. A similar divergence is 
possible while the subpopulations maintain a limited spatial contact [3 [T71 HH1 HHJ 120] , 
which is commonly referred to as parapatric speciation. Pre- or post-mating incompati- 
bilities can develop subsequently, leading to reproductive isolation. Neither the primary 
phenotypic separation nor the secondary reproductive isolation require disruptive selec- 
tion. Hence, allopatric or classical parapatric speciation may well take too much space 
and time to account for radiation events and rapid species divergence |2"T] I22j. 

In recent years, phylogeographic observations have produced convincing evidence for 
speciation in sympatry. Reproductive isolation has occurred in cichlid populations in 
African lakes over a few thousand generations PU 121 0] ■ A salmon population is reported 
to have separated within only 14 generations [22] • Sympatric speciation thus appears to 
occur rapidly even in small contiguous environments without spatial barriers. In theoret- 
ical models, it is always driven by disruptive selection. A phenotypic split can be favored, 
for example, if individuals of similar phenotype compete more strongly than distant ones 
(HI E] ■ In a sexually reproducing population, however, such splits can only happen if the 
subpopulations become reproductively isolated so that the birth of hybrids is suppressed 
|3 13 El HHj • Of course, the sympatric scenario cannot explain the spatial population struc- 
ture observed in the phylogeographic studies. Spatial patching of subpopulations appears 
to be ubiquitous. For example, the sister species of cichlids tend to organize themselves 
into neighboring regions [31 HI 121] . Another well-documented case are phytophagous in- 
sects, which are found to evolve mating assortativity together with specificity to different 
host plants [2*5] . 

These observations call for a more comprehensive model that captures the divergence in 
phenotypic traits, mating, and space as a cooperative dynamical process. Only recently 



IBM simulations in an extended model space with environmental fitness gradient have 
been presented extending previous studies with complete spatial mixing [T51 ITi] . 

The model discussed here addresses parapatric speciation, i.e., generic intermediate 
cases between sympatry and allopatry. It affords a detailed analysis of the dynamics, 
allowing us to identify different primary speciation mechanisms and their conditions of 
occurence. (A concise discussion of classification issues and of the appropriate terminology 
can be found in Ref. [2T].) It turns out that the basic evolutionary forces driving speciation 
can be captured by a deterministic "react ion- diffusion" approach. We also discuss the role 
of stochastic effects as they appear in individual-based models. In this way, we recover 
the well-known mechanisms of allopatric and sympatric speciation. However, there are 
many environments with inhomogeneities on smaller scales in space and time (such as 
in the examples quoted above), where spatial variations prevent sympatry and diffusive 
migration prevents allopatry. Adaptative evolution then operates by a new mechanism, 
for which we propose the term diapatric speciation. The population reaches a final state 
of efficient spatial patchiness and phenotypic different at ion without hybrids, which is 
triggered and sustained by assortative mating. This is in contrast to the traditional 
view of parapatric speciation, where assortative mating takes a merely secondary role in 
reinforcing an existing boundary between emerging species |17| H%1 I21j. 

Model 

We consider a population via its density in "internal" and "external" space, iV = iV(x; r; t). 
Internal coordinates x = (xi, . . . ,x n ) denote phenotypic quantities, e.g., body size, beak 
length, colour. Internal coordinates can be inherited. This representation is purely phe- 
notypic. A comparison with explicitly genetic models is given below. 

External coordinates r = (?"!,..., r^) lie in the simplest case in <i-dimensional Eu- 
clidean space. More complicated geometries, e.g., network structures of habitat patches 
in fragmented landscapes, are also possible. In this work, we focus on habitats with a 
gradient in quality for different phenotypes, which induce a spatial dependence of the 
optimal phenotype x opt (r) and a population density iV(x, r; t) with a joint dependence on 
internal coordinates x and external coordinates r. 




Figure 1: The fitness landscape of a heterogeneous model environment involves a fitness funtion 
f(x, r) that depends on a trait variable x and a spatial coordinate r. (a) The left region (r < 0) 
favors smaller values of x, the right region (r > 0) larger ones. The optimal trait x op t(r) varies 
between the values ±xo over a spatial interval given by the total size r^. (b) At a given point 
r, the fitness is maximal at x op t(r) and decays rapidly over a characteristic scale Wf, called the 
niche width. 

In the simplest version of the model we consider one phenotypic coordinate and a 
one-dimensional external space of size tq, so N = N(x, r; t) whith — tq/2 < r < tq/2. The 



phenotype x is directly related to an ecological fitness or carrying capacity, e.g., with the 
explicit choice 

/(*, r) = /(* - * opt (r)) = exp (- {% - (1) 

which is taken to be constant in time. It decreases with the distance of x from x opt = 
xq sin(7rr/r ), on a scale Wf in phenotype space, xq is a measure for habitat heterogeneity 
and r is the spatial scale of variation. For an illustration see Figure ^ 
The population N(x, r; t) is subject to the dynamics 

d t N(x, r; t) = Xd 2 N(x, r; t) + R(x, r; t) + (f(x, r) - K(x, r; t)) N(x, r; t) (2) 

which has the form of a react ion- diffusion equation. 

The simplest type of motion in the population is diffusion, in Eq. (J2J) appearing as the 
term Xd 2 N, to which we restrict ourselves in this work. The prefactor defines a length 
scale in space, r\ = vA, which has to be compared with the habitat size r . 

The special case of eq. (J2J) with R = describes the dynamics of an asexual or clonal 
population. It similar to the familiar Lotka-Volterra form. The resource supply f(x,r) 
and the competition load 

K(x, r;t) = J dy (3{x, y) N{y, r; t), (3) 

which sums up the influence of individuals of trait y on those with trait x, combine to 
the frequency-dependent fitness f — K. The competition kernel 

P(x, y) = (3(\x - y\) = exp Ul^lli) (4) 

is maximal for x = y and decays on a scale wp in internal space. 

Extending this approach to sexually reproducing populations requires a more detailed 
model for birth processes, whose rate itself becomes dependent on the maternal and 
paternal population densities. It is convenient to introduce the birth excess per phenotype, 
space, and time 

R(x,r;t) = dy dz C(x\y, z) m(y, z;t) N(z,r;t) — N(x,r;t). (5) 



by summing over the density of possible mothers N(z, r; t) multiplied by the probability 
density m(y, z) for a z-female to mate with a y-male and the inheritance probability den- 
sity C(x\y,z) that this couple will have offspring of phenotype x. The subtracted term 
N(r,t) describes the total birth rate in the clonal limit. With the definitions of C and m 
given below, it is easy to check that / dx R(x, r) = 0. Hence, the excess birth rate describes 
the net reshuffling of population density through sexual reproduction, and / — K remains 
a useful measure of the frequency- and space-dependent fitness. The genetic function C 
is approximated by a Gaussian, C(x\y,z) = exp(— (x — x) 2 / (2wc{x) 2 )) / \j2irwc(x) 2 , with 
x = (y+z)/2, so offspring is distributed near the mean of the parents' phenotype. More- 
over the standard deviation wc{x) changes little over the relevant range of pheno types. 
This form can be justified from the hypergeometric model (23122112011301, where the quan- 
titative trait x is encoded by L independent two-allele loci with equal allele frequencies. 
However, provided the number of independent loci is sufficiently large, it remains valid 
more generally, even if (i) the number of loci changes or (ii) the symmetry between the 
loci is lost jUI] because allele frequencies change or linkage disequilibria delevop during 



the speciation process. Typically this would result in a decrease of Wc, but as long as 
wc < Wf and wq < w@, variations in wq do not influence the results significantly. See 
also the discussion at the end of this Section where we show that this form of C(x\y,z) 
emerges from a genetically explicit model quite generically. 

Mating preference is crucial for the development of any structure in the population. 
Without it, the mating probability is just proportional to the available males. In this case 
the entire population is mixing and forms a single cluster in phenotype, see Figure a). 
This changes with an affinity of females towards certain types of males, 

/ .x v(y,z)N(y,r;t) 

Here we restrict ourselves to assortative mate choice by the ecological trait within a range 
of width 

(\x — y\^ \ 
2 • ( 7 ) 

With strong enough mating assortativity reproductively isolated subpopulations can co- 
exist, as shown in Figure E^b). 

The population dynamics (J2J) always leads to a stationary density N(x,r), which 
reflects the primary selection given by the fitness functions / and K. On longer, evolu- 
tionary time scales, the population evolves through secondary selection, i.e., by adaptive 
mutations modifying its mating range We study this process starting from a 

spatially uniform initial state with random mating. A single step involves an initially 
small mutant population that invades the resident population and eventually becomes a 
new stationary state N(x,r) with different trait and mating characteristics. At each step 
we evaluate whether a stationary state N(x,r) with given is unstable with respect to 
a small mutant population n(x,r;t) with different mating range. Successful mutants are 
found to invade the resident population completely, producing a new stationary state. A 
possible dependence w^(x) due to a linkage disequilibrium, not taken into account here, 
is expected only to enhance the selection pressure towards assortativity. If adaptive sub- 
stitutions are sufficiently rare, an evolutionary pathway can be represented as a sequence 
of intermediate stationary states leading to an evolutionary stable final state N es (x,r) 
[55] . Along the pathway, the number of adaptive steps parametrizes evolutionary time. 
More generally, the mating range may be thought of as a further quantitative trait, 
the population state being described by a joint distribution N(x, w^, r). The distribution 
of is strongly peaked, which justifies the approximation of Eq. (J2J). The average value 
of W/j, evolves along fitness gradients towards the final state. Generic evolutionary stable 
states are found to have either random or strongly assortative mating. 

Interesting variations in the internal structure of the model are related to the mating 
preference. It can depend on ecologically neutral but inheritable traits such as mating 
time, marker traits, and in all cases one may observe phenotypic differentiation ^21 El 
I n smaller populations some individuals may be unable to mate. Assortativity 
restricts the number of possible mates and should be disfavored under such circumstances. 
With some modification of Eq. (0) this effect can be studied and it turns out that certain 
types of reproductive isolation are actually favored [TOJ EH] . 

Unlike in our model inheritance in sexual population dynamics is often modelled ge- 
netically explicit. The classical approach is to consider a locus with two alleles, say a 
and A, and under which conditions preferentially homozygous subpopulations develop 
0. In computer simulations longer "genomes" can be used, typically two strings of L 



bits with the "alleles" and 1. Genome space is then very large, 2 , and a common 
way to follow the evolution of a population are simulation of so called individual based 
models (IBMs), [T3*l ITU l2*o]. For their evaluation population characteristics are sampled 
over large populations, long times, and many independent runs. 

Based on phenotypes but closely related to genetics is the so called hypergeometric 
model O I2E1 123 IHO] , where the phenotype of an individual with 2L loci is a quantitative 
trait proportional to the number of one type of alleles, e.g. 

2L 

x = ^^G{0,l,...,2L}, (8) 

and the alleles are o v G {0, 1}. If all genotypes mapping onto a phenotype are equally 
probable in a population, one can derive the probability C(x\yz) for a couple with phe- 
notypes y and z to have offspring with x: explicitly for a haploid and to a very good 
approximation for a diploid genome [2E1 12Z|- Going one step further away from the un- 
derlying genetic concept leads to models of Quantitative Genetics j21| one of which is 
ours. 

Generally such models neglect gene fixation. Also the hypergeometric model |T2J 123 
123 EH! ma y be invalidated as the central assumption of equiprobability of the various 
genotypes contributing to one phenotype can fail |3*T| . But the same difficulty also arises 
for IBMs, as e.g. in [T3J H3J |2H| , where only a "good" choice of mutation rate, population 
and genome size allows for meaningful dynamics with respect to the question of speciation. 
It is in these cases, that the phenotype related hypergeometric model and also quantitative 
phenotypic models as ours behave similarly and thus remain meaningful. 

In Figure |2]we show some examples of C(x\y, z) as functions of the offspring's pheno- 
type x for fixed phenotypes of mother and father, obtained by sampling over an IBM with 
genome length 2L = 64, population size 16384, random mating, run for 10 7 generations. 
Phenotypes are given by Eq. (JHJ). The children's phenotypes are distributed around values 
x(y, z) = (y + z)/2 with a (nearly Gaussian) distribution whose widths wc are practically 
independent of the parents' phenotypes. In panel (a) are examples for three values of 
z = y, in panel (b) for values z = 2L — y, such that the parents' mean phenotypes are 
all identical x — L. Panel (c) shows a long time average of the proportion of 1-alleles in 
the entire model genome compared to the values observed at the 2L model loci. These 
simulations show that the elementary combinatorial rules of inheritance on the genome 
level used in typical 1MB simulations can quite well be approximated on the phenotype 
level by continuous functions for C(x\y, z). The maximum value x and the width wq may 
be subject to corrective terms, but the principle structure of C(x\y,z) remains valid. 

In fact it turns out that the precise functional form of the interactions does not mat- 
ter too much. We have also studied alternative forms of faster or weaker decay, e.g., 
f(x,r) ~ exp(— (x/wf) 4 ), (3(x) ~ exp(— (x/wp) 2 ), //(x) ~ exp(— x/w^). Important are 
the length scales in internal and external space: the inheritance uncertainty Wc, the 
competition range wp, the resource width Wf, habitat heterogeneity xq, extent of habi- 
tatvariation ro, and migration range r\. Their combination and mutual relation decides 
about reproductive separation of the population into two or more subpopulations. 

Results 

We first discuss the special case where f(x,r) = f(x) and N(x,r;t) = N(x;t) do not 
depend on the spatial coordinate r, which requires Xq = 0, and we consider a spatially 
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Figure 2: Examples of C(x\y, z) sampled from a simulation of 16384 randomly mating individ- 
uals with genome length 2L = 64 over 10 7 generations with mutation probability 10~ 3 per locus 
and generation, (a) Examples for y = z = 20, 32, and 44. b) Same z, but y = 2L — z. (c) Allele 
frequencies of l's at the single loci differ by less than 1% from the average value 1/2. 



averaged "mean field" version of the model. Also in this limit sympatric speciation can 
become manifest, analogous to the results of the IBM in Refs. |1H[ IT4"] which is set in 
a similar ecological frame. For all evaluations we assumed the ecological interactions to 
extend over a wider range than the inheritance uncertainty, wc < Wf and Wc < wp. 
During its adaptation the assortativity range u> M varies, but it remains larger than wc- 

Evolving reproductive isolation can lead to separation into subpopulations. Equilib- 
rium profiles of Eq. (J2J) are shown in Fig. EJ in panel (a) a unimodal population for random 
mating, in panel (b) a bimodal under mating assortativity after evolutionary adaptation 
of w^. For a large enough relative width of the habitat Wf/wp > 1.1 the population 
evolves into a speciating state as in (b), otherwise it remains unseparated. 

For given parameters %, wp, Wf there may be a range of w^, where both a unimodal 
(as in Fig. Ota)) and a bimodal population (as in Fig. Efb)) are stable fixed points of 
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Figure 3: Fixed point configurations for Equation (J2J). (a) weak mating preference gives uni- 
modal population structure, (b) strong mating preference allows for bimodal population struc- 
ture. The phenotypic scale is indicated by the resource curve f(x) plotted in light gray, the 
vertical scale of population density is arbitrary. 

Eq. (J2J). As assortativity gets stronger (decreasing w^) the unimodal profile ( Fig. 0(a)) 
becomes less stable. It is interesting to note that the transition between a unimodal 
and a bimodal population density is not a gradual process but a fast switch at a critical 
assortativity range w^. The switch occurs on the time scale of population dynamics, much 
faster than the evolutionary adaptation of w^. Under conditions where the evolutionary 
stable can increase again (e.g. slow variation of Wf) one finds hysteresis between the 
jumps from uni- to bimodal populations and back. 

A population profile as in Fig. 0(b) cannot be a stable fixed point of the asexual version, 
i.e., the limit — > 0, of Eq. (J2J): The gap between the two parts of the population would 
fill up resulting in a wider unimodal population profile covering most of the accessible 
phenotype range ^Hj. When u> M remains finite two peaks in the population profile having 
widths close to but mutual distance greater than are stabilized by sexual reproduction 
because it accumulates offspring closer to their maxima. By this effect sexual reproduction 
helps speciation. 

We now turn to the general case. The spatial model dynamics generates differentiation 
of the population, which can be measured in two ways: (i) The mating differentiation index 
5 is defined as the actual rate of crossmating between two subpopulations at a given point 
r, normalized by the same rate with random mating. For any population state N(x,r), 
the mating differentiation index 5 (at the point r = and between the subpopulations 
x < and x > 0) is defined by 

1-6 = 1= [ dyf dzl[m(y,z,0)N(y,0) + m(z,y,0)N(z,0)], (9) 

Zj Jy<0 Jz>0 I 

where Z is the same integral evaluated with random mating, i.e., with fi(y,z) = 1 for 
all y, z. (Analogous measures can be defined for different r and different trait subpop- 
ulations). Here we monitor the two subpopulations x < and x > at the boundary 
between the left and right regions (r = 0). (ii) The spatial differentiation index x is 
defined in terms of the "trait overlap" between the populations at two different points in 
space. Phenotypes x that are intrinsically viable at one of these points are distinguished 
from those that are merely advected by diffusion. The spatial differentiation index \ 
(evaluated at the points — r /2 and ro/2) is defined by 

1 - X = ~r, 1 /oW , 1 ... fdxN v (x,-r /2)N v (x,r /2) (10) 
Z(r /2) Z(-r /2) J 

where Z(r) = fdx N v (x, r). A phenotype x is counted as intrinsically viable at the point r 
if wc/3(0)N(x, r) 2 — Xd 2 N(x, r)/dr 2 > 0, i.e., if a small nonzero population in the interval 



[x — wc(x)/2,x + wc(x)/2] could exist even without diffusive advection. In this case, 
we set N v (x,r) = N(x,r), otherwise N v (x,r) = 0. Here we take the points r = — r /2 
in the left region and r = tq/2 in the right region. Both indices vary between (no 
differentiation) and 1 (complete separation). 

Following the differentiation in phenotype and space over evolutionary times, three 
main mechanisms can be identified. They are distinguished by the structure of their 
evolutionary stable final populations N es (x,r), measured, for example, by the resulting 
differentiation indices 5 es and Xes- 
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Figure 4: Three mechanisms of speciation can be distinguished by the time dependence of the 
mating differentiation index 5 and the spatial differentiation index % (see text). The initial 
population has random mating (5 = 0) and is spatially homogeneous (x = 0). Primary selection 
with random mating (left part of the diagrams) is followed by secondary selection on the mating 
range w m (right part of the diagrams), (a) Allopatric speciation: Continuous evolution by 
primary selection towards spatial segregation (x es = 1) without reproductive isolation (8 es < 1). 
(b) Sympatric speciation: Discontinuous evolution towards reproductive isolation (<5 es = 1) 
without spatial segregation (x es < 1). (c) Diapatric speciation: Discontinuous, cooperative 
evolution towards reproductive isolation (5 es = 1) and spatial segregation (x e s = !)■ 
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Figure 5: Evolutionary stable populations after speciation. Left column: Populated regions in 
the (x,r) plane, given by N cs (x,r) > 0. Intrinsically viable phenotypes (shown in black) are 
distinguished from those advected by diffusion (grey). Right three columns: trait distributions 
N es (x,r = — tq/2) (left region), N es (x,r = 0) (boundary between left and right region), and 
N es (x,r = ro/2) (right region), (a) Allopatric speciation: One contiguous population cluster, 
unimodal trait distributions, species boundary with hybrids, (b) Sympatric speciation: Two 
disjoint clusters, bimodal trait distributions, (c) Diapatric speciation: Two disjoint clusters, 
trait distributions unimodal within the regions and bimodal at the boundary, species boundary 
without hybrids. 
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Figure 6: Phase diagram of speciation, specifying the mechanism as a function of the effective 
environment heterogeneity xo/wf, the effective competition range wp/wf, and the diffusion 
coupling between the regions, A/Vq. (a) Cross-section in the variables wp/wf and xq/wj at 
fixed A/Vq = 0.01. (b) Cross-section in the variables A/Vq and xq/wj at fixed wp/wj = 2.0. 
Diapatric speciation is the generic mechanism in heterogeneous environments with diffusive 
coupling. 

Allopatric speciation shows a gradual increase of the spatial differentiation index \ 
up to Xes = 1, see Fig. Ufa). This expresses the patching of small-x phenotypes into the 
left region and large-x phenotypes into the right one. The spatial adaptation of traits 



involves primary selection by the fitness function f(x, r) only and occurs independently of 
mating behavior. Since there is no sufficient selection pressure towards assort at ivity, the 
mating differentiation S es remains small. The corresponding evolutionary stable popula- 
tion N es (x,r) is a contiguous cluster in the (x,r) plane as shown in Fig. E(a). At a given 
point r, the trait distribution is unimodal and centered around the local fitness maximum 
Xf(r). There is a limited gene flow between the large- x and small- x subpopulations, which 
is maintained by the intermediate phenotypes near the boundary (r = 0). 

Sympatric speciation is characterized by an increase of mating differentiation up to 
complete reproductive isolation; see Fig. ^b). The index S jumps discontinuously from a 
value So < 1 to S es = 1, implying that stationary population states with So < S < 1 cannot 
exist. The speciation is driven by secondary selection involving the frequency-dependent 
fitness K, just as in previous models of strict sympatry. Spatial variations are irrelevant, 
and the spatial differentiation remains incomplete (xes < 1)- The evolutionary stable 
population N es (x,r) shown in Fig. Efb) consists of two disjoint clusters, corresponding to 
a bimodal trait distribution at every r. The gene flow between these subpopulations is 
suppressed by assortative mating. In particular, there are no hybrids near the boundary 
(r = 0). 

Diapatric speciation is the co-evolution of assortative mating and spatial segregation 
by secondary selection. The indices S and x jump to S es = Xes = 1 simultaneously, 
leading to an evolutionary stable state with reproductive isolation and patching into the 
left and right region, see Fig.HJc). Prior to the jump, the spatial segregation is prevented 
by diffusive coupling between the regions. It becomes possible only once reproductive 
isolation is established. The population N es (x, r) of Fig. EJc) has two disjoint clusters. 
The trait distribution is unimodal within both regions and bimodal near the boundary; 
there are again no hybrids. The suppression of the gene flow between the clusters is now 
two-fold, by reproductive isolation and by spatial separation. 

Does a population actually speciate, and if so, by which mechanism? This turns out to 
depend largely on only three parameters, the effective environmental heterogeneity Xo/wf, 
the effective competition range Wp/wf, and the diffusive coupling between the regions, 
X/tq. (Here we have chosen wj and r as the basic scales in trait space and real space.) 
The "phase diagram" of Fig. El shows the mechanism of speciation as a function of these 
parameters. Allopatric speciation is possible only with a sufficiently large heterogeneity 
and a sufficiently small diffusive coupling (i.e., small A or large region size tq). Conversely, 
sympatric speciation requires a sufficiently small heterogeneity, as well as a sufficiently 
small competition range (wp/wj < 1). Diapatric speciation involves no restriction on 
the competition range, that is, it works for frequency-dependent as well as for density- 
dependent selection. It is seen to be the generic mechanism in many realistic environments 
with heterogeneities and diffusion. 

This compares to the recent results of Ref. where adaptive speciation is seen to 
generate a sharp geographical differentiation in an individual based model. Working over 
a wider range of parameters the present model is able to relate this diapatric mechanism 
to other modes of speciation by identifying their respective regions in terms of the relevant 
parameters. 

The present model thus allows for a clear identification of the evolutionary mecha- 
nisms underlying speciation, of the dynamical patterns, and of the resulting population 
structures. Examples are the separation indices S and x (Eqs. JHJ) and (JTOJ) ) and the dis- 
tinction of intrinsically viable populations from populations merely advected by diffusion. 
Of course, this kind of differential analysis is very difficult in individual-based models, 
which always suffer from small discrete population sizes. On the other hand, the effect 



of demographic stochasticity and other fluctuations can also be studied within the frame- 
work of Eq. by adding a stochastic noise term. The evolutionary stable population 
densities iV es (x, r) are found to be stable under such perturbations. Stability or instability 
of stationary states N(x, r) become immediately apparent in differential equations such 
as by their rates of convergence or divergence. The fast transition between coherent 
and segregated population states thus explains itself naturally from a simple analysis. 

Discussion 

In summary, our model suggests that speciation is a highly cooperative process involving 
the adaptive differentiation of a population in its ecological characters, its mating be- 
havior, and its spatial structure. Diapatric speciation is the generic mechanism of fully 
coupled differentiation. It reduces to allopatric or sympatric speciation in special cases. 
All three mechanisms are part of a unified dynamical picture, conceptually different from 
the old dichotomy between sympatry and allopatry. 

In classical observations, the spatial separation of newly formed species has often been 
regarded as the primary driving force of the speciation process. Our results call for a 
fresh look at the data and may offer a different interpretation in some cases. The diapatric 
mechanism involves spatial separation as an induced process, triggered by the reproductive 
isolation of subpopulations. This two-fold separation in phenotype and space between the 
emerging species cuts the gene flow more efficiently than the other mechanisms, which 
involve only one kind of separation. 

Diapatric species boundaries are established and maintained by natural selection so 
no external barriers have to be postulated. They follow regional boundaries and are dis- 
tinguished from the allopatric case by the efficient suppression of hybrids in the boundary 
zone during the primary speciation process. (Secondary reproductive isolation can sup- 
press hybrids also in allopatry.) It is crucial to note that reproductively decoupled popu- 
lations can adapt to spatial heterogeneities of smaller size than interbreeding ones. Diap- 
atric speciation is also remarkably fast, since the loss of interbreeding takes place through 
an abrupt change of the stationary population state as discussed above. This transition 
is driven by natural selection, unlike the secondary mating differentiation mechanisms 
in allopatry, which are expected to operate by genetic drift and hence to be slower. Of 
course, the genetic fixation of permanent incompatibilities between the emerging species 
(postzygotic isolation) is always slow. Before that point, both reproductive and the spatial 
separation are reversible if the environmental conditions change, as has been confirmed 
by recent observations jSH]- Hence, diapatric splits appear to be an efficient adaptation 
mechanism for sexual populations on small scales of space and time. Most of these splits 
are wiped out again on longer time scales, while a few develop into permanent speciation. 
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